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Abstract 

The effect of the Kerr nonlinearity on linear non-diffractive Bessel beams is inves- 
tigated analytically and numerically using the nonlinear Schrodinger equation. The 
nonlinearity is shown to primarily affect the central parts of the Bessel beam, giving 
rise to radial compression or decompression depending on whether the nonlinear- 
ity is focusing or defocusing, respectively. The dynamical properties of Gaussian- 
truncated Bessel beams are also analysed in the presence of a Kerr nonlinearity. It 
is found that although a condition for width balance in the root-mean-square sense 
exists, the beam profile becomes strongly deformed during propagation and may 
exhibit the phenomena of global and partial collapse. 
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1 Introduction 



Diffractive spreading of waves is a classical phenomenon in wave dynamics and 
an inherent feature of beam propagation. Much attention has been devoted 
to the possibility of counteracting the dispersive spreading by focusing effects 
due to medium nonlinearities e.g. the Kerr effect, cf. [1]. However, it has 
also been pointed out, [2], that non-diffracting beams are possible also in 
linear media. In particular, the Helmholz equation that governs the linear 
diffractive dynamics of a wave beam allows classes of diffraction-free solutions. 
In addition to the plane wave solutions, the two-dimensional counterparts, the 
cylindrically symmetric Bessel solutions, also propagate with preserved form, 
while also allowing for a concentrated beam profile. The drawback from an 
application point of view is the fact that these beams have infinite energy, 
and consequently cannot be realized physically. Various ways to circumvent 
this problem have been suggested, the most obvious being to truncate the 
Bessel beam at some radius e.g. by a Gaussian truncation, forming the so 



Preprint submitted to Elsevier Science 



2 February 2008 



called Bessel-Gauss beams, [3]. While such a truncation clearly reintroduces 
diffraction, the beam broadening could be made small if the propagation length 
is kept smaller than the corresponding diffraction length of the Bessel-Gauss 
beam. In particular, since the Bessel beam diffracts sequentially, starting with 
the outer lobes, cf . [4] , the central part of the beam remains intact for a certain 
distance of propagation. 

Recently, there has been growing interest in nonlinear effects in connection 
with Bessel and Bessel-Gauss beams, [5,6,7,8]. Of special interest for the 
present investigation is the attention given to media with an intensity de- 
pendent refractive index, i.e., Kerr media, see e.g. [7]. The work carried out 
in [7] considers the limit of weak nonlinearity, which makes it possible to use 
a perturbation approach involving an expansion around the lowest order lin- 
ear (stationary) Bessel solution for solving the evolution equation, being the 
nonlinear Schrodinger equation. 

In the present paper, we investigate in more detail the nonlinear generalisa- 
tion of the linear diffraction-less Bessel beam solutions as well as the nonlinear 
dynamics of Bessel-Gauss beams. Stationary solutions in two dimensions are 
determined by the Bessel equation modified by a nonlinear term, i.e., the radi- 
ally symmetric nonlinear Schrodinger equation. The modified Bessel solutions, 
the "nonlinear Bessel beams", are studied using approximate analytical and 
numerical methods. The results show that the nonlinearity primarily affects 
the central high intensity parts of the beam profile, which become radially 
compressed or decompressed depending on whether the nonlinearity is fo- 
cusing or defocusing, respectively. The beam profile for large radii remains 
a Bessel function with a phase shift being the only remaining effect of the 
nonlinearity. However, for the defocusing nonlinearity an amplitude threshold 
exists, above which no solutions decaying to zero exist. 

The dynamical properties of Gaussian-truncated Bessel beams in the presence 
of a Kerr nonlinearity are also studied. An exact analytical solution was pre- 
viously found for the linear dynamics of the Bessel-Gauss beams, [3]. Based 
on the virial theorem, which gives an exact analytical description of the vari- 
ation of the beam width in the root-mean-square (RMS) sense, important 
information about the effect of the nonlinearity on the beam dynamics can 
be obtained. In particular, it is found that a focusing nonlinearity tends to 
cause an evolution stage where the central parts of the Bessel beams are 
initially compressed. Depending on the strength of the nonlinearity, different 
scenarios are possible, e.g. the subsequent evolution may involve an essentially 
diffraction-dominated behaviour, but for increasing nonlinearity, two forms of 
collapse may appear. Either a part of the beam collapses, while the RMS width 
of the beam still increases, or above a certain threshold, the RMS width goes 
to zero in a finite distance. Numerical simulations of the dynamics illustrate 
the different scenarios. 
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2 The nonlinear Schrodinger equation 



The propagation of an optical wave in a nonlinear Kerr medium is described 
by the nonlinear Schrodinger equation. This implies that the slowly varying 
wave envelope, ip(z, r), of a cylindrically symmetric beam satisfies the following 
equation 



where z is the distance of propagation, fc is the wave number in vacuum, and 
k is the nonlinear parameter. Additional physical effects, e.g. attenuation and 
gain, can be modelled by using complex coefficients in Eq. (1). The obtained 
complex equation, which in one dimension has analytical soliton solutions, [9], 
is the cylindrical generalisation of the Pereira-Stenflo equation. It has been 
investigated using a variational approach, [10], but further work is needed 
to fully describe the complex case. In the present work, the coefficients are 
assumed to be real. It is convenient to introduce the normalisation f = r/a , 
where a is a characteristic width of the beam, and z = z/L D , where L D = 
2/c ao is the Rayleigh length. Eq. (1) then takes the form 



.dip = d 2 ip | Idj) | . 2 
dz df 2 f df 



where k = L d k. For simplicity we will suppress the tilde in the subsequent 
expressions. 

We begin the analysis by looking for stationary solutions of Eq. (2). For this 
purpose we write ip = tp(z, r) = A(r) e lSz , which leads to the eigenvalue equa- 
tion 

(PA IdA r A ,o 

— + -— + 5A + kA 3 = 0. 3 

ar z r ar 



This equation is to be solved subject to the boundary conditions that the so- 
lution should be finite when r = and go to zero as r — > oo. The lowest order 
solution in the physical situation when the nonlinearity balances the diffrac- 
tion, i.e., the focusing case with k > 0, corresponds to the so called Townes 
soliton [11], which has essentially the same properties and sec/i-shaped form as 
the one-dimensional soliton solution, cf. [12]. In particular, this solution only 
exists for negative eigenvalues, which are uniquely related to the maximum 
amplitude. 



3 



The Bessel beams are the solutions of the linear Schrodinger equation (k = 0) 
and are given by A(r) = A Jq(\/5 r). Clearly, well-behaved solutions exist 
only for positive eigenvalues 5. In contrast to the nonlinear case, the linear 
eigenvalue problem has a continuous set of (positive) eigenvalues, which are 
independent of the amplitude of the beam profile. The first task of the present 
analysis is to analyse the nonlinear Bessel beams, being the solutions of Eq. (3) 
for positive 5 and k ^ 0. We note that by introducing 

\5\r, A=y/W8\A, (4) 



only the signs of 5 and k remain in Eq. (3). Thus, without loss of generality, 
it will be assumed that 5=1 and k — ±1. 



3 Nonlinear Bessel beams 



In order to analyse the properties of the nonlinear Bessel beams analytically 
it is instructive to start by examining the central part of the pulse, which is 
determined by \J~5r <C 1. Since the initial derivative of A(r) must be zero, 
we have the approximation A(r) = A + 0(r 2 ), which implies that to second 
order in r, Eq. (3) can be approximated by the linear equation 

^ + if + (S + KA})A ^, (5) 



with the corresponding solution 

A = A J (^5 + KA 2 r). (6) 

This solution is valid for small r only, but nevertheless gives important infor- 
mation about the nonlinear modifications of the Bessel beam. In the focusing 
case (k > 0), the main lobe tends to be compressed and we expect that the 
amplitude of the nonlinear Bessel beam will oscillate more rapidly than the 
linear Bessel beam. On the other hand, in the defocusing case (k < 0), the 
main lobe should become wider and the nonlinear solution should oscillate 
slower than the linear one. In particular, if the amplitude is chosen to fulfil 

5 + kA 2 < 0, (7) 



the expression under the square root will be negative. This corresponds to 
the modified Bessel function, which is growing with r. Thus the presence of a 
defocusing nonlinearity can qualitatively change the behaviour of the solution. 
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It is clear that if the second derivative of A is positive initially, it will remain 
positive. This can be seen by rewriting the equation as 




(8) 



where 8 e g = 8 + kA 2 . If <5 c fr is negative at r = 0, the solution will be growing 
for small r and 5 c s will be further decreased. This increases the derivative of 
the solution, and implies that a sufficiently strong defocusing nonlinear term 
will give rise to a monotonically growing solution, which is not compatible 
with the condition at infinity. 

A more accurate description of the main lobe of the solution can be obtained 
using variational analysis and Ritz optimisation, cf. [12]. When using varia- 
tional analysis, it is important to find a good set of trial functions that gives 
tractable calculations while maintaining sufficient accuracy. A trial function 
that should approximate the main lobe of the nonlinear Bessel beam reason- 
ably well is 



where jo is the first zero of the Bessel function. This also has the advantage 
that the exact linear result is recovered in the limit k — > 0. In the variational 
procedure, we assume r to be given and consider A as a free parameter. The 
Lagrangian corresponding to Eq. (3) is 




(9) 




(10) 



o 



where 





Using the ansatz (9) we obtain 



£ = \ [io^iU) - rl5AlJl{j G ) - C.rlKAi 




where 



i 




(13) 



o 
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The variation with respect to A , i.e., dC/dA = yields 




(Jo 2 



-r 2 6)Jf(j ) 





(14) 



5J?(j ) + 2C lK A 2 



By setting k = 0, the linear result \/5r = j is recovered. If k is positive the 
value of Tq is decreased, which corresponds to compression of the main lobe. 
A negative k gives the opposite effect. This result confirms the previously 
obtained picture of the effects of the nonlinearity on the main lobe of the 
linear Bessel solution. 

The result of the variational analysis is compared with numerical solutions in 
Figs. 1 and 2. Different r -values have been used, and they are easily identi- 
fied in the figures, since the variational approximation is zero when r = tq. 
For clarity the plotted curves have been normalised with respect to their am- 
plitudes at r = 0. It is seen that the Bessel ansatz represents a rather good 
approximation in the focusing case, but that the presence of the nonlinearity 
changes the shape for small r, making it more peaked than the Bessel pro- 
file. In the defocusing case, the nonlinear solution instead has a flatter form 
than the linear Bessel function, Fig. 2. It is also seen that for increasing r , or 
equivalently increasing A , the approximation deteriorates. This is due to the 
fact that in this case there exists a threshold value for the amplitude in order 
to have well-behaved solutions, cf. (7). The critical value for the amplitude 
is Aq = 1. The variational result also predicts this behaviour, although, as is 
inferred from Eq. (14), the critical value is found to be slightly different 



Clearly, as A(0) approaches the threshold value, which is unity, we expect the 
accuracy of the variational approximation to deteriorate. 

We now turn to an investigation of the overall behaviour of the nonlinear 
Bessel beam profiles. The fact that the main influence of the nonlinear term 
is a rescaling of the radial coordinate makes it reasonable to look for an ap- 
proximate solution of the form 



where f(r) is a function of r, 5, k, and A . It is difficult to determine / using 
analytical methods, but by noticing that the linear solution can be written as 




(15) 



A(r) » A J {f(r)), 



(16) 




(17) 
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Fig. 2. Comparison of the variational results (solid lines) to the numerical ones 
(dashed lines) for a defocusing nonlinearity. 

and by comparing with Eq. (5) it seems reasonable that a good approximation 
should be obtained by the implicit expression 

A = A J (J y/5 + KA 2 (r')dr^j . (18) 
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Although this is, in fact, an integral equation for A(r), it nevertheless provides 
a very simple formula for finding A numerically. The corresponding approx- 
imate solution is compared to the numerical solution of the full equation in 
Fig. 3. When the amplitude is low the two curves are identical, since the ansatz 
then reduces to the Bessel function. In the case of a focusing nonlinearity, there 
is good agreement between the two approaches, but it is also seen that a phase 
shift appears between the curves for increasing Aq. Quite good agreement is 
seen also in the defocusing case. In particular, the initial flattening is well 
modelled. The phase shift is now of the opposite sign. 

This approximate solution implies that the argument of the Bessel function 
increases approximately as J r yS + kA 2 dr', which is a nonlinear generalisation 
of the linear case. Thus the main effect of a focusing nonlinearity is to increase 
the curvature of the peaks by increasing the growth rate of the argument, 
making the solution radially compressed. In the defocusing case the curvature 
is decreased, which is most clearly seen in the main lobe. 
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Fig. 3. The implicit analytical solution given by Eq. (18) (solid lines) together with 
the numerical result (dashed lines). The different initial amplitudes are indicated in 
the graph. A defocusing nonlinearity is used in the fourth plot. 

Finally, Figs. 4 and 5 further illustrate the nonlinear deformations of the linear 
diffraction-less Bessel solutions. The numerically obtained curves clearly show 
the features discussed above; the radial compression of the central lobe in the 
focusing case and the radial expansion in the defocusing case. The expansion 
effect in the latter case rapidly increases as the amplitude approaches the 
critical value A — 1, above which no stationary solutions are possible. The 
phase shifting effect of the nonlinearity on the Bessel-like oscillations is also 
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seen, the shift changing sign with the sign of the nonlinearity. 




Fig. 4. A focusing nonlinear term gives rise to a radial compression, which is illus- 
trated using numerical simulations. 
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4 Analysis of nonlinear Bessel-Gauss beams 



The linear diffraction properties of Bessel-Gauss beams have been analysed 
and solved analytically, see [3]. In the present section we will analyse the 
nonlinear dynamics of beams, which initially have a profile in the form of 
a Bessel function truncated by a Gaussian. Since a general solution of this 
problem cannot be given, we will use the virial theorem to obtain analytical 
information and numerical simulations for determining the evolution of the 
beam profile. 

The virial theorem, see e.g. [13,14,15], provides exact and explicit information 
about the dynamic variation of the width of the beam, a, defined in the RMS 
sense as a 2 = (r 2 ), where 



</(r)> = 



^f{r)\jj{z,r)\ 2 rdr 
JS°\1>(z,r)\*rdr ' 



(19) 



The virial theorem asserts that 



d 2 a 2 
dz 2 



H 



constant, 



(20) 



where / and H are invariants of the two-dimensional nonlinear Schrodinger 
equation, Eq. (1), and are defined as follows 



J = y \ip(z, r)\ 2 r dr, 



H 



dip(z, r) 



dr 



r dr. 



(21) 
(22) 



This means the invariants correspond to the (integrated) beam intensity and 
the Hamiltonian. Thus, the virial theorem implies that a 2 must be a second 
order polynomial in z, with coefficients determined by the initial beam profile, 
■0(0, r). For initial phase functions that do not depend on r, the linear term 
in z vanishes, and the beam width varies as 



(z)=a 2 (0) ^l + sign^)-^ 



(23) 



10 



where L is a characteristic length given by 

2 



J' 2 



AH 



cr 2 (0)7 



4^ 



dip(z,r) 
dr 



lz,r) 



r dr 



/ °°r 2 |^(0,r)| 2 r dr 



(24) 



Clearly this approach cannot be used for analysing the linear, or the nonlin- 
early modified, stationary Bessel beam solutions of the nonlinear Schrodinger 
equation since all integrals involved in the virial theorem are infinite. However, 
for a physical beam, with finite integral content, the virial theorem is useful. 
In general it is seen that with weak nonlinear focusing effects, the Hamiltonian 
is positive and the RMS width will increase quadratically with a characteristic 
diffraction length given by L . When the amplitude of the beam increases, the 
Hamiltonian decreases and eventually changes sign. This implies that the RMS 
width goes to zero after a finite length equal to Lq — the well known nonlinear 
collapse phenomenon, where L now plays the role of collapse length. 



For Bessel-Gauss beams, [3] , the initial profile is of the form 



V>(0,r) = A J (^j exp ^~ 



(25) 



Inserting this into Eq. (24) we obtain the following expression for the charac- 
teristic length L 



r 
2L ( 



sm-asm 

Ss(p) 



(26) 



where p = r\l p\, A = kA^Tq, and the integrals S n , n = 1,2, 3, are given by 



00 

Si(p) = J [pxJ (x) + J\{x)} 2 exp(— /j,x 2 )x dx, 


00 

S2(p) = J Jo ( x ) exp(— 2fix 2 )x dx, 


00 

S , 3 (/i) = J x 2 Jq(x) exp(—fix 2 )xdx. 



(27) 
(28) 
(29) 



Since we are primarily interested in the case when the Gauss function truncates 
the outer parts of the Bessel function, we have p = r 2 jp\ <^ 1. In this limit, 
the asymptotic values of the integrals, S n , n = 1,2, 3, are obtained analytically 

as 
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3 

S 2 (p)^D 1 - -j-j In//, 



4tt 2 
1 



4v^F// 3 / 2 ' 



(30) 

(31) 
(32) 



Since S 2 (p) goes rather slowly towards infinity as p becomes small, it is nec- 
essary to determine the constant D\ in order to have good accuracy for finite 
p. Using numerical evaluation of the integral we find D\ 0.202. This implies 
that the characteristic length can be approximated as 



' r 

,2Lr 



2p 



(33) 



with D 2 ~ 0.715. In the linear case, the characteristic length is seen to scale 
simply as Lq oc po, i.e., the diffraction is determined solely by the truncation 
radius and as p 00 , the non-diffracting Bessel beam is recovered. For 
increasing values of the nonlinearity parameter, A, but for a fixed truncation 
radius, the value of L increases and for a certain critical value of A, the 
nonlinearity balances the diffraction to give a Bessel-Gauss beam, which is 
diffraction-less in the RMS sense. 



5 Dynamics of Bessel-Gauss beams 



When a truncated linear Bessel-Gauss beam propagates, the diffraction ini- 
tially affects only the outermost part of the pulse, where the truncation has 
changed it from the Bessel shape. The central parts are initially diffraction bal- 
anced and remains so until the "diffraction front" propagating inwards from 
the outer parts eventually reach the inner lobes and also these parts start 
to diffract outwards. On the other hand, the nonlinear effect is strongest at 
the centre of the beam, where the intensity is highest, and with a focusing 
nonlinear term, the main lobe will start to compress. In fact, it will start to 
compress irrespective of the degree of nonlinearity since the linear diffraction 
is already balanced. If the nonlinear effect is weak the compression will even- 
tually stop, and diffraction will become the dominating effect. This evolution 
is illustrated in Fig. 6, where an FDTD simulation using A — 1 and p = 0.01 
is shown. Although the central parts initially compress, the virial theorem pre- 
dicts beam broadening in the RMS sense. Clearly this is no contradiction since 
the broadening of the outer parts more than compensate the compression of 
the centre. In order to further illustrate how the main lobe is compressed, the 
intensity at r = has been plotted for different initial amplitudes as a func- 
tion of propagation distance in Fig. 7. The curves have been normalised with 
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respect to their amplitudes at r = and 11 = 0.01. We emphasise the oscillat- 
ing behaviour for the highest amplitudes. If the nonlinear effect is sufficiently 
strong, the virial theorem predicts a collapse of the beam to zero RMS width. 
According to Eq. (33), the amplitude threshold for this type of behaviour is 
A > A c ps 2.26 for ii = 0.01. It is well known that the nonlinear evolution of 
two-dimensional beams may lead to a break up of the beam into a diffracting 
background profile with a monotonously compressing filament, that collapses 
in a finite distance of propagation. Thus, the width of the filament goes to zero 
and the intensity becomes infinite whereas the beam width in the RMS sense 
still increases. The simulations for the present case of Bessel-Gauss beams 
show that when the amplitude is increased further above Aq = 1, the small 
second peak of Fig. 7 will start to dominate and eventually the simulations 
indicate a collapse of the second peak, although the RMS width still tends to 
infinity. In fact, even if the RMS width remains constant, the beam should 
still be able to undergo partial collapse. 

Much effort has been devoted to the study of two-dimensional collapse phe- 
nomena induced by the Kerr nonlinearity, see e.g. [13,14,15] and references 
therein. In particular, it has been found that the virial theorem poses a suf- 
ficient but not necessary condition for the occurrence of a singularity where 
the amplitude becomes infinite. Thus the appearance of a partial collapse sin- 
gularity below the threshold for a global collapse, as predicted by the virial 
theorem, is in accordance with earlier results. 
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Fig. 7. The intensity at r = as a function of the propagation distance using 
different initial amplitudes. 

6 Conclusions 

Based on the nonlinear Schrddinger equation in cylindrical geometry we have 
studied the modification of the diffraction-less linear Bessel beams caused by 
the nonlinear Kerr effect. The stationary as well as the dynamic properties 
of the solutions to this equation have been analysed and both analytical and 
numerical techniques have been used. The investigation shows that the non- 
linearity primarily affects the main and inner lobes of the Bessel beams. In 
the case of the stationary solutions, the central region of the nonlinear Bessel 
beam tends to become radially more narrow or more extended depending on 
whether the nonlinearity is focusing or defocusing, respectively. Asymptoti- 
cally the solutions are still of the same oscillating form as the diffraction-less 
Bessel beams, the only remaining feature of the nonlinearity being a phase 
shift as compared to the linear case. However, in the case of the defocusing 
nonlinearity, there is a finite amplitude threshold for well-behaved solutions to 
exist. Above this limit the nonlinear diffraction effect becomes larger than the 
linear effect, which counteracts the diffraction, and no solutions are possible 
which vanish at infinity. 

The properties of Gaussian-truncated Bessel beams have also been studied in 
the presence of the Kerr nonlinearity. It has been shown, using the virial the- 
orem, that a non-diffracting situation in the RMS sense is possible to obtain 
by balancing nonlinear focusing and linear diffraction. However, this situation 
does not correspond to a stationary case of the beam profile. Significant re- 
distribution of the beam occurs and using numerical simulations, the dynamic 
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interaction between linear diffraction and nonlinear focusing has been anal- 
ysed for varying degrees of nonlinearity. It has been found that, in particular 
the central parts of the beam may become significantly distorted and may even 
partially collapse even though the beam width, defined in the RMS sense, re- 
mains constant or even increases. This result is in agreement with the classical 
picture of dynamic self focusing of two-dimensional beams in nonlinear Kerr 
media. 
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